Toolkit

# In order to get ggmap running, run the following:
install.packages("devtools")
library(devtools)
devtools::install_github("dkahle/ggmap")
devtools::install_github("hadley/ggplot2")
library(tidyverse) # all of our usual tools
library(ggmap)     # for spatial viz.
library(viridis)   # for color scales

Oklahoma earthquakes

ok_quakes <- read.csv("OKquakes.csv")
str(ok_quakes)
'data.frame':   5204 obs. of  9 variables:
 $ time     : Factor w/ 5204 levels "2000-10-08T10:16:23.780Z",..: 5204 5203 5202 5201 5200 5199 5198 5197 5196 5195 ...
 $ latitude : num  35.8 36.6 35.7 36.8 35.7 ...
 $ longitude: num  -97.6 -98.8 -97.4 -97.8 -97.4 ...
 $ depth    : num  6.48 6.29 5.97 5.41 6.18 ...
 $ mag      : num  2.6 3.2 2.5 3 2.8 2.5 3 3.2 3.1 3.3 ...
 $ place    : Factor w/ 1433 levels "0km E of Pawnee, Oklahoma",..: 466 840 1177 1069 1236 318 1236 1092 1242 90 ...
 $ year     : int  2015 2015 2015 2015 2015 2015 2015 2015 2015 2015 ...
 $ month    : int  12 12 12 12 12 12 12 12 12 12 ...
 $ type     : Factor w/ 4 levels "light","minor",..: 4 2 4 2 4 4 2 2 2 2 ...

Displaying spatial data with ggplot2

ggplot(data = filter(ok_quakes, year == 2015), aes(x = longitude, y = latitude)) +
  geom_point()

The ggmap package

# Download the requested map
ok_map <- get_map(location = 'oklahoma', zoom = 7)
Map from URL : http://maps.googleapis.com/maps/api/staticmap?center=oklahoma&zoom=7&size=640x640&scale=2&maptype=terrain&language=en-EN&sensor=false
Information from URL : http://maps.googleapis.com/maps/api/geocode/json?address=oklahoma&sensor=false
# Display the map
ggmap(ok_map)

  1. Look at the help file for get_map and read about the location and zoom arguments. Do the following:

There are other options aside from a terrain map.

  1. Look at the help file for the get_map function and read about the maptype and source arguments. Specify different maptypes and sources. How do things change?

Now that we have our map of Oklahoma, we can add points to it:

ggmap(ok_map) +
  geom_point(data = ok_quakes, aes(x = longitude, y = latitude, color = year))

A better idea might be to facet the plot by year. In order to to this we need to add the base_layer argument to the ggmap function:

ggmap(ok_map, base_layer = ggplot(ok_quakes, aes(x = longitude, y = latitude))) +
  geom_point(size = .8) +
  facet_wrap(~year, ncol = 5)

  1. Add color to the faceted map based on the type of earthquake. Be sure to use an appropriate color palette. Do any patterns emerge?

  2. The file CAquakes.csv contains data on earthquakes in California from January 1, 2010 through the end of 2015. Use ggmap to plot the earthquakes that occurred in California during this time, similar to what we did for Oklahoma above. Do any patterns emerge?

Heatmaps

ok2015 <- filter(ok_quakes, year == 2015)

We can create heat maps using stat_density2d in place of our geom_point layer:

ggmap(ok_map, base_layer = ggplot(ok2015, aes(x = longitude, y = latitude))) +
  stat_density2d(aes(fill = ..level.., alpha = ..level..), size = 0.01, 
    bins = 16, geom = "polygon") +
  scale_fill_viridis() +
  scale_alpha(guide = FALSE)

  1. Create a heat map of California earthquakes in 2015.
LS0tCnRpdGxlOiAiQ3JlYXRpbmcgbWFwcyBpbiBSIgphdXRob3I6ICJDTVNDIDIwNSwgV2ludGVyIDIwMTciCmRhdGU6ICJNYXJjaCA1LCAyMDE3IgpvdXRwdXQ6IAogIGh0bWxfbm90ZWJvb2s6CiAgICBjc3M6IGxhYi5jc3MKLS0tCgoKIyMgVG9vbGtpdAoKYGBge3IgZXZhbD1GQUxTRX0KIyBJbiBvcmRlciB0byBnZXQgZ2dtYXAgcnVubmluZywgcnVuIHRoZSBmb2xsb3dpbmc6Cmluc3RhbGwucGFja2FnZXMoImRldnRvb2xzIikKbGlicmFyeShkZXZ0b29scykKZGV2dG9vbHM6Omluc3RhbGxfZ2l0aHViKCJka2FobGUvZ2dtYXAiKQpkZXZ0b29sczo6aW5zdGFsbF9naXRodWIoImhhZGxleS9nZ3Bsb3QyIikKYGBgCgoKYGBge3IgbWVzc2FnZT1GQUxTRX0KbGlicmFyeSh0aWR5dmVyc2UpICMgYWxsIG9mIG91ciB1c3VhbCB0b29scwpsaWJyYXJ5KGdnbWFwKSAgICAgIyBmb3Igc3BhdGlhbCB2aXouCmxpYnJhcnkodmlyaWRpcykgICAjIGZvciBjb2xvciBzY2FsZXMKYGBgCgojIyBPa2xhaG9tYSBlYXJ0aHF1YWtlcwoKLSBUaGUgbnVtYmVyIG9mIGVhcnRocXVha2VzIGluIE9rbGFob21hIGhhcyBpbmNyZWFzZWQgaW4gdGhlIHBhc3QgZGVjYWRlCgotICBXZSB3aWxsIGV4YW1pbmUgZGF0YSBvYnRhaW5lZCBmcm9tIHRoZSBVLlMuIEdlb2xvZ2ljYWwgU3VydmV5IChVU0dTKSBvbiBlYXJ0aHF1YWtlcyBpbiBPa2xhaG9tYSBzaW5jZSAyMDAwIAoKCmBgYHtyfQpva19xdWFrZXMgPC0gcmVhZC5jc3YoIk9LcXVha2VzLmNzdiIpCnN0cihva19xdWFrZXMpCmBgYAoKIyMgRGlzcGxheWluZyBzcGF0aWFsIGRhdGEgd2l0aCBgZ2dwbG90MmAgey5zbWFsbGVyfQoKYGBge3J9CmdncGxvdChkYXRhID0gZmlsdGVyKG9rX3F1YWtlcywgeWVhciA9PSAyMDE1KSwgYWVzKHggPSBsb25naXR1ZGUsIHkgPSBsYXRpdHVkZSkpICsKICBnZW9tX3BvaW50KCkKYGBgCgoKIyMgVGhlIGBnZ21hcGAgcGFja2FnZQoKLSBDb250ZXh0IGlzIGltcG9ydGFudCwgd2UgbmVlZCBtb3JlIHRoYW4ganVzdCB0aGUgcG9pbnRzIQoKLSBgZ2dtYXBgIGFsbG93cyB1cyB0byBncmFiIHRoaXMgY29udGV4dCBmcm9tIEdvb2dsZSBvciBvdGhlciBzb3VyY2VzCgotIEJ5IGRlZmF1bHQsIGBnZXRfbWFwYCwgZG93bmxvYWRzIGEgdGVycmFpbiBpbWFnZSBmcm9tIEdvb2dsZSBtYXBzCgpgYGB7cn0KIyBEb3dubG9hZCB0aGUgcmVxdWVzdGVkIG1hcApva19tYXAgPC0gZ2V0X21hcChsb2NhdGlvbiA9ICdva2xhaG9tYScsIHpvb20gPSA3KQoKIyBEaXNwbGF5IHRoZSBtYXAKZ2dtYXAob2tfbWFwKQpgYGAKCjEuIExvb2sgYXQgdGhlIGhlbHAgZmlsZSBmb3IgYGdldF9tYXBgIGFuZCByZWFkIGFib3V0IHRoZSBgbG9jYXRpb25gIGFuZCBgem9vbWAgYXJndW1lbnRzLiBEbyB0aGUgZm9sbG93aW5nOgotIFVzZSBgbG9jYXRpb24gPSAnbGF3cmVuY2UgdW5pdmVyc2l0eSdgIHRvIGdldCBhIG1hcCBvZiB0aGUgc3Vycm91bmRpbmcgYXJlYS4gVHJ5IGRpZmZlcmVuY2UgdmFsdWVzIGZvciBgem9vbWAuCi0gQ2hvb3NlIGEgY2l0eSBhbmQgdXNlIFtHb29nbGUgbWFwc10oaHR0cHM6Ly93d3cuZ29vZ2xlLmNvbS9tYXBzKSB0byBmaW5kIHRoZSBsb25naXR1ZGUgYW5kIGxhdGl0dWRlIG9mIGEgbGFuZG1hcmsuIFNwZWNpZnkgdGhlIGNvb3JkaW5hdGVzIGluIHRoZSBgbG9jYXRpb25gIGFyZ3VtZW50IGFuZCBwcm9kdWNlIGEgbWFwIG9mIHRoZSBzdXJyb3VuZGluZyBhcmVhLiBBZ2FpbiwgdHJ5IGRpZmZlcmVudCB2YWx1ZXMgZm9yIGB6b29tYC4KClRoZXJlIGFyZSBvdGhlciBvcHRpb25zIGFzaWRlIGZyb20gYSB0ZXJyYWluIG1hcC4KCgoyLiAgTG9vayBhdCB0aGUgaGVscCBmaWxlIGZvciB0aGUgYGdldF9tYXBgIGZ1bmN0aW9uIGFuZCByZWFkIGFib3V0IHRoZSBgbWFwdHlwZWAgYW5kIGBzb3VyY2VgIGFyZ3VtZW50cy4gU3BlY2lmeSBkaWZmZXJlbnQgYG1hcHR5cGVzYCBhbmQgYHNvdXJjZXNgLiBIb3cgZG8gdGhpbmdzIGNoYW5nZT8KCk5vdyB0aGF0IHdlIGhhdmUgb3VyIG1hcCBvZiBPa2xhaG9tYSwgd2UgY2FuIGFkZCBwb2ludHMgdG8gaXQ6CgpgYGB7cn0KZ2dtYXAob2tfbWFwKSArCiAgZ2VvbV9wb2ludChkYXRhID0gb2tfcXVha2VzLCBhZXMoeCA9IGxvbmdpdHVkZSwgeSA9IGxhdGl0dWRlLCBjb2xvciA9IHllYXIpKQpgYGAKCkEgYmV0dGVyIGlkZWEgbWlnaHQgYmUgdG8gZmFjZXQgdGhlIHBsb3QgYnkgeWVhci4gSW4gb3JkZXIgdG8gdG8gdGhpcyB3ZSBuZWVkIHRvIGFkZCB0aGUgYGJhc2VfbGF5ZXJgIGFyZ3VtZW50IHRvIHRoZSBgZ2dtYXBgIGZ1bmN0aW9uOgoKYGBge3J9CmdnbWFwKG9rX21hcCwgYmFzZV9sYXllciA9IGdncGxvdChva19xdWFrZXMsIGFlcyh4ID0gbG9uZ2l0dWRlLCB5ID0gbGF0aXR1ZGUpKSkgKwogIGdlb21fcG9pbnQoc2l6ZSA9IC44KSArCiAgZmFjZXRfd3JhcCh+eWVhciwgbmNvbCA9IDUpCmBgYAoKMy4gQWRkIGNvbG9yIHRvIHRoZSBmYWNldGVkIG1hcCBiYXNlZCBvbiB0aGUgYHR5cGVgIG9mIGVhcnRocXVha2UuIEJlIHN1cmUgdG8gdXNlIGFuIGFwcHJvcHJpYXRlIGNvbG9yIHBhbGV0dGUuIERvIGFueSBwYXR0ZXJucyBlbWVyZ2U/Cgo0LiAgVGhlIGZpbGUgYENBcXVha2VzLmNzdmAgY29udGFpbnMgZGF0YSBvbiBlYXJ0aHF1YWtlcyBpbiBDYWxpZm9ybmlhIGZyb20gSmFudWFyeSAxLCAyMDEwIHRocm91Z2ggdGhlIGVuZCBvZiAyMDE1LiBVc2UgYGdnbWFwYCB0byBwbG90IHRoZSBlYXJ0aHF1YWtlcyB0aGF0IG9jY3VycmVkIGluIENhbGlmb3JuaWEgZHVyaW5nIHRoaXMgdGltZSwgc2ltaWxhciB0byB3aGF0IHdlIGRpZCBmb3IgT2tsYWhvbWEgYWJvdmUuIERvIGFueSBwYXR0ZXJucyBlbWVyZ2U/IAoKCiMjIEhlYXRtYXBzCgpgYGB7cn0Kb2syMDE1IDwtIGZpbHRlcihva19xdWFrZXMsIHllYXIgPT0gMjAxNSkKYGBgCgoKV2UgY2FuIGNyZWF0ZSBoZWF0IG1hcHMgdXNpbmcgYHN0YXRfZGVuc2l0eTJkYCBpbiBwbGFjZSBvZiBvdXIgYGdlb21fcG9pbnRgIGxheWVyOgoKYGBge3J9CmdnbWFwKG9rX21hcCwgYmFzZV9sYXllciA9IGdncGxvdChvazIwMTUsIGFlcyh4ID0gbG9uZ2l0dWRlLCB5ID0gbGF0aXR1ZGUpKSkgKwogIHN0YXRfZGVuc2l0eTJkKGFlcyhmaWxsID0gLi5sZXZlbC4uLCBhbHBoYSA9IC4ubGV2ZWwuLiksIHNpemUgPSAwLjAxLCAKICAgIGJpbnMgPSAxNiwgZ2VvbSA9ICJwb2x5Z29uIikgKwogIHNjYWxlX2ZpbGxfdmlyaWRpcygpICsKICBzY2FsZV9hbHBoYShndWlkZSA9IEZBTFNFKQpgYGAKCjUuIENyZWF0ZSBhIGhlYXQgbWFwIG9mIENhbGlmb3JuaWEgZWFydGhxdWFrZXMgaW4gMjAxNS4=